Role of inertia in two-dimensional deformation and breakup of a droplet 
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We investigate by Lattice Boltzmann methods the effect of inertia on the deformation and break- 
up of a two-dimensional fluid droplet surrounded by fluid of equal viscosity (in a confined geometry) 
whose shear rate is increased very slowly. We give evidence that in two dimensions inertia is necessary 
for break-up, so that at zero Reynolds number the droplet deforms indefinitely without breaking. 
We identify two different routes to breakup via two-lobed and three-lobed structures respectively, 
and give evidence for a sharp transition between these routes as parameters are varied. 



The role of inertia in the deformation and breakup 
of a droplet in a shear flow remains an open problem 
in fluid mechanics. In three dimensions, a droplet in a 
steady shear field will break at a critical capillary num- 
ber Ca ~ 0.5 even without inertia 0. Crudely speaking, 
this happens when it is deformed enough to undergo the 
Rayleigh instability, which is the peristaltic instability of 
a long tube of one static fluid in another, driven by in- 
terfacial tension. Ca measures the relative importance of 
viscous to interfacial forces; the effect of inertia, whose 
importance is quantified by the Reynolds number Re, is 
essentially perturbative. Thus, although situations arise 
where a droplet that would be stable at zero Re is un- 
stable at a finite Re, this really amounts to an inertia- 
induced shift of the critical Ca. The contribution of iner- 
tia was studied in recent numerical work ^ [| ||, |f| . Much 
of this work uses a flow protocol based on sudden onset 
of shearing at different flow rates, rather than a gradual 
ramping up of the rate; this saves much time numerically 
but describes a distinct physical situation from the slow 
ramp studied here; both are of interest experimentally. 

The situation is somewhat different in two dimensions. 
Although not realizable directly in laboratory experi- 
ments, this case is more accessible by simulation || 
and is interesting because the nonlinear physics involved 
in breakup may not be the same as in 3D. In 2D, the 
Rayleigh instability is effectively switched off, because 
periodic variations in the width of the droplet (at fixed 
volume) always increases the amount of interface. In the 
present work we study numerically in 2D the breakup of 
a droplet, and find quite distinctive nonlinear physics, 
in which inertia plays an essential part. The relatively 
moderate demands of 2D simulations allows us to estab- 
lish the scenario relatively clearly, even for the case of an 
'infinitesimal' ramp, which we here approach numerically 
by a novel recursive updating scheme for the flow rate. 

Method: We use the lattice Boltzmann method to sim- 
ulate a (symmetric) binary fluid deep within the two- 
phase region. The phase separation of the order parame- 
ter <f> is induced by a chemical potential p — A(j> + B<f? — 
kV 2 </> derived from a <j) 4 Landau free energy The 
evolution equation for the order parameter is 



with M a mobility; this is coupled with a Navier-Stokes 
equation 



p(d t u + u. Vu) = -Vp + 4>Vfx + ryV 2 u 



(2) 



where u is the fluid velocity, p its density and r\ is the vis- 
cosity. The interfacial tension is included in this model 
through the gradient terms in the chemical potential and 
it is given by a 



MV — 2kA Note that topological re- 
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connections of the fluid-fluid interface are handled implic- 
itly by the order parameter diffusion and no singularity is 
encountered at breakup. Care was taken with parameter 
selection to ensure that the droplet behaviour is domi- 
nated by hydrodynamic and capillary forces, except very 
close to the breakup point, where both diffusion, and an 
interaction arising from the finite breadth of the inter- 
faces (a few lattice spacings), kick in. 

In our simulations, we first set up a droplet of fluid in 
a surrounding immiscible phase of the same viscosity. To 
impose a shear flow with strain rate 7, we introduce mov- 
ing walls at the top and bottom edges of the simulation 
cell with opposite velocities ±Lj/2 where L is the cell 
width. A droplet in such a shear profile will deform with 
the flow until and unless the viscous stress is balanced by 
interfacial tension. The ratio of these forces is expressed 
in the capillary number 



Ca = 



(3) 



where Rq is the radius of the undcformcd droplet; the 
Reynolds number measures the ratio of inertial and vis- 
cous forces: 



Re 



(4) 



We quantify the deformation of a droplet by calculat- 
ing the following moments of the order parameter field 



Or 



u.V(/> = MV fi 



(1) 



b = y,°(w ( 5 ) 

X 

c a = ^9((f))x a (6) 

X 

d a p = y]6(<f>)<f>(x a - c a ){xf3 - cp) (7) 
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Here b gives the 2D droplet 's area and c/b the position 
vector of its centre. The eigenvalues of d a p/b give us two 
length scales l\ > I2 corresponding to the two axes of the 
droplet. The Taylor deformation is denned as 



Note that very deformed droplets have D close to (but 
still less than) unity. 

To study the breakup behavior of droplets and find 
the critical capillary number (beyond which no steady 
state for a single droplet stably exists) we need to find 
the stationary shapes of droplets at various Ca and Re. 
We ramp up the shear rate in a way that finds all stable 
droplet shapes before the breakup. To do this efficiently 
we use an algorithm that increases the shear rate 7 step- 
wise and then waits for the droplet shape to reach steady- 
state before calculating Ca and Re. If there is breakup, 
or if the resulting deformation increment is too large (ex- 
ceeding a threshold AD C ~ 0.05) a smaller increase A7 
in the shear rate is tried. This process is iterated until 
either equilibrium is achieved, or A7 falls below a pre- 
set limit chosen very close to zero (Ajc ~ 10~ 7 ). In the 
latter case, the last value of 7 for which equilibrium was 
reached is identified as the final stable shear rate; this 
fixes the value of Ca and Re at breakup for the chosen 
run. By varying simulation parameters (such as viscos- 
ity) an ensemble of curves of Re vs Ca was generated. 
Each termininated in a breakup point. 

The simulations reported below have lattices bounded 
by moving walls separated by 120 lattice sites, with pe- 
riodic boundary conditions along the flow direction (in 
which the lattice length is 400). Other parameters com- 
mon to all runs were, in lattice units —A = B = 
0.06; k = 0.039, from which the interfacial tension is de- 
rived as 0.046; and M — 3.0. Ten viscosity values were 
selected in the range 0.23 — 0.88, again in lattice units ||. 
These parameter values were chosen to minimise lattice 
anisotropy, spurious diffusive currents, and other well- 
documented but fairly controllable artefacts of the lattice 
Boltzmann method j|, [To). The initial droplet radius in 
each run was 21 lattice units. However, it was found 
that the increased curvature of a sheared droplet could 
lead to slight changes in its size as a result of increased 
Laplace pressure causing a shift in solubility. This effect 
was minimized by using the true size of the droplet when 
calculating Rq for the purposes of Re and Ca (see below). 

The ratio L/Rq ~ 5.7 is not so large as to approach 
the behaviour of a droplet being sheared in an infinite 
medium. Wall proximity effects are present in our results 
(as they are in many experiments), but we do not study 
these systematically here. In addition, the slight change 
in Rq with shear rate mentioned above means that the 
wall effects vary slightly during each run. By varying 
lattice size and shape we checked that there was no undue 
influence on the results. The same applies to the effect 



FIG. 1: For zero Reynolds number a two dimensional droplet 
does not break up. It extends until the width of the droplet 
is of the order of the interface width. The picture shows a 
droplet with Ca = 1.1; the upper wall is moving rightward, 
the lower, leftward. 



caused by proximity of a droplet to its periodic image 
along the flow direction. 

Although it is hard to quantify precisely all systematic 
errors in lattice Boltzmann Q we believe the quantitative 
results reported below are generally good to within ten 
percent and mostly rather better. Our main conclusions, 
however, are qualitative in nature. 

Zero Reynolds number results: First we report a 
simulation for zero Reynolds number. While lattice 
Boltzmann fully treats inertial effects, a stationary zero 
Reynolds number flow can also be simulated (essentially 
by switching off the term in the algorithm that handles 
fluid convection). In contrast to previous suggestions 
H [HJ, we found that a two-dimensional droplet in a 
shear flow at zero Re does not break up by any fluid- 
mechanical process. Instead it continues to deform (Fig- 
ure [l]) until the width of the droplet is comparable to 
the interface width, at which point dissolution-assisted 
breakup occurs. This corresponds to a molecular rather 
than a hydrodynamic breakup mechanism. A similarly 
elongated droplet in 3D would certainly not survive, due 
to the Rayleigh instability. 

Finite Reynolds number results: Both the Reynolds 
number and the capillary number depend linearly on the 
shear rate. We can, however, form the dimensionless 
number / = Re/Ca which will be almost constant for 
each run (deviations arising only from the small changes 
in droplet size Rq; see above). The nominal / values 
for the nine different runs simulated are identified in the 
caption to Figure ||, where we show deformation curves 
in each case. At modest Ca values the curves with a 
higher / (more inertia) show a larger deformation. A 
likely explanation for this is the Bernoulli effect. This 
will create regions of lower pressure at the edges of the 
droplet, tending to stretch it into the flow so that it expe- 
riences a stronger flow than applies for a droplet at zero 
Re. 

The simulation corresponding to the lowest curve is the 
zero Reynolds number case and here the droplet does not 
break, as previously discussed. All other curves end at 
the last stable configuration for the droplet that could 
be found using our ramp method. It can be difficult to 
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FIG. 2: Taylor deformation vs. capillary number under 
slow ramping of shear rates for (top to bottom) nine runs 
of viscosities in the range 0.181 — 0.877 (lattice units) and 
the zero Re limit. Values of the parameter I — Re/Ca are 
I = 18, 13, 11, 9.1, 7.3, 6.3, 3.7, 2.2, 1.3, (top to bottom). 

establish the stability of the very slowly evolving, near- 
critical droplets that arise close to this last stable configu- 
ration (particularly if they are very deformed) so we have 
taken care to run the simulations for tens of thousands 
of timesteps per iteration in these parts of the curves. 

With increasing / we find that there is a trend towards 
steeper deformation curves, and at a critical I c between 
6.3 and 7.3 we observe a singularity in which the defor- 
mation curve appears to develop an infinite slope at its 
end point. Further deformation curves with larger / all 
end similarly at a deformation of about 0.6. We per- 
formed simulations for a number of different lattice sizes 
and the exact value of I c depends on the cell dimensions, 
with increasing wall separation causing a shift to smaller 

As explained earlier, our ramp algorithm is designed 
to allow only a certain amount of deformation for each 
new stationary state. So, before the algorithm reaches a 
region of infinite slope it will force increasingly smaller 
shear-rate increments; when these fall below the pre-set 
threshold the algorithm will save the last configuration 
and exit. To check the fate of these droplets we reloaded 
the last configurations manually and again increased the 
shear rate by a very small amount. All these droplets 
then broke up into two droplets; none into three. 

We show the dependence of the last stable capillary 
number on Re in Figure |3|. Each of the runs appears an 
(almost) straight line on the plot, ending at the breakup 
point. Two different critical capillary number regimes 
are clearly separated in the Re/Ca plane. Figure || gives 
examples of the last stable droplet shapes for both I < I c 
and I > I c . For I < I c we find that the droplets develop 
a double dimple; beyond the critical capillary number 
Ca c (7) these break up into three. For I > I c the last 
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FIG. 3: Deformation curves plotted in the Reynolds number / 
capillary number plane. Note that the zero Re curve coincides 
with the x-axis and does not end in breakup. The envelope 
of breakup points is shown to guide the eye: this has two 
disconnected parts, corresponding to two-lobe and three-lobe 
droplets at the last stable configuration. Typical shapes are 
as indicated. 

stable droplet shows one dimple; beyond Ca c these break 
up into only two droplets. 

Note that, in these simulations, we only observed 
breakup for Re < 1.2. The shape of the breakup en- 
velope in Figure ||, together with the knowledge that 
the zero Reynolds number case does not lead to break- 
up, suggests that the envelope could have a horizontal 
asymptote at some nonzero value of Re (perhaps close 
to unity). This would imply a finite critical Re below 
which no breakup could occur, however large Ca. How- 
ever, more detailed work would be needed to reach a firm 
conclusion on this point. 

We also performed simulations with a less sophisti- 
cated ramp algorithm in which the shear rate was in- 
creased as before, but without restriction on the defor- 
mation increment. In this case, different behaviour was 
seen around I c . Droplets which, when ramped more 
carefully, would show the two-lobed breakup mode (with 
a vertical tangent to the deformation curve at its end 
point; Figure |^) could now jump past this onto the up- 
per branch of the deformation curve, developing instead 
a stable three-lobed shape, finally breaking at higher Ca 
into three droplets. The deformation curve for one such 
case is compared to the more careful ramp for the same 
/ in Figure ^. This proves that the details of the flow 
history can strongly affect the final breakup mode and 
also the shear rate at which it occurs. 

Comparing this result with the general trend shown 
for larger / in Figure suggests a conjecture: that the 
behaviour is governed by an underlying, continuous de- 
formation curve that has, for I > I c , folded over into an 
S-shape. With such a curve, one might naturally expect 
a discontinuous jump in deformation on ramping up to 
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FIG. 4: Two deformation curves for I — 11. For one run 
with the deformation restriction was removed. The droplet 
overshoots the underlying break-up point and jumps into a 
three lobe structure, which survives to larger capillary num- 
bers until eventually breaking into three droplets. 



to a three lobed structure before breaking, no matter how 
slow the ramp rate, and then breaks into three droplets. 
From extrapolation of the critical capillary number we 
proposed that, even at very large Ca, a threshold in Re 
of order unity may need to be exceeded; but an alterna- 
tive is that the Ca c diverges only as Re-^ 0. 

Because of the care required in achieving equilibration 
close to Ca c it would be difficult to perform similar nu- 
merical studies in three dimensions without very substan- 
tial resources. Nonetheless, an important open question 
is whether, on increasing the relative importance of iner- 
tia, the breakup mode in three dimensions also undergoes 
a sudden transition as we discovered here. We are un- 
aware of experiments addressing this, and note that in 
such experiments (as well as in simulations) very careful 
control of the shear history might be needed to identify 
any transitions that may be present. 
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the point of vertical tangency. However, it appears that 
droplets cannot survive the jump (but instead break in 
two) if brought very close to this point before the jump 
is made. With a faster ramp, the jump occurs sooner; 
the droplet then can survive the jump and continue onto 
the upper branch of the deformation curve. The mecha- 
nism seems to involve transient stretching just after the 
shear rate is incremented; this allows the droplet to find 
the upper branch of the deformation curve, but only if 
the shear rate increment is sufficient. Note that such an 
S-shaped deformation curve would also admit hysteresis, 
in which two states of deformation were both possible at 
exactly the same Ca. However, we never got clear evi- 
dence of this in our simulations, so any hysteresis loop 
appears to be relatively small. 

Conclusions: In this article we have provided evi- 
dence that droplet breakup in two dimensions requires a 
nonzero Reynolds number. The dependence of the crit- 
ical capillary number Ca c on the dimensionless number 
7= Rc/Ca for two-dimensional droplet break-up shows 
an unexpected transition at I — I c . (The precise value 
of I c depends on the wall separation.) This transition 
is consistent with an evolution of the underlying defor- 
mation curve _D(Ca) into an S-shape. (However we did 
not observe hysteresis around this transition.) A slow 
enough ramp leads to a two-lobed droplet breakup mode 
at Ca c for I > I c whereas under different ramp condi- 
tions a droplet can, before reaching Ca c , jump from a 
two lobed structure to a three lobed one, which is then 
stable to much higher Ca. For I < I c the droplet evolves 
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